Monocyte analysis V6

---

Analysis of monocyte data (SCZ vs CTR): simplest model (~ status)

Load data

First step QC

Plotting the first QC results:

1. Density plot

2. Covariate correlation

#Covariate correlation
for (i in colnames(int.covs))
  int.covs[[i]] <- as.numeric(int.covs[[i]])
MEs_lin(int.covs,int.covs)

3. PC-covariate correlation

#PC-covariate correlation (nomal p)
PCA_cov_cor_R(int.covs,trnsf.df)

Next step: QC on vst-normalized expression values.

1. Gene expression variance (low variance in lowly expressed genes (to the left) = good)

meanSdPlot(assay(vsd))

2. PCA

#Check influence of covariates on data variance after transformation
PCAplot(pca, "batch", PoV.df=PoV, colors=c('#2E7281','black','#50B8CF',"green"), pc.1 = 3, pc.2 = 4)
PCAplot(pca, "status", PoV.df=PoV, pc.1 = 3, pc.2 = 4)

DESeq2 anaylsis

We then go into the differential gene-expression analysis:

1. Overview of the results

out of 24437 with nonzero total read count adjusted p-value < 0.1 LFC > 0 (up) : 26, 0.11% LFC < 0 (down) : 15, 0.061% outliers [1] : 0, 0% low counts [2] : 0, 0% (mean count < 0) [1] see ‘cooksCutoff’ argument of ?results [2] see ‘independentFiltering’ argument of ?results

As well as the number of differentially expressed genes at lfc </> -0.5/0.5 at padj < 0.1:

[1] 36

Interactive results table

1. Results

#DT::renderDT(data.frame(res.complex), "OH1.5A",scrollY=1000)

createDT(data.frame(res), "Status", scrollY=1000)

Plotting results

2. Volcano plot:

EnhancedVolcano(res,
                lab = rownames(res),
                x = 'log2FoldChange',
                y = 'padj',
                pCutoff = 0.1,
                FCcutoff = 2,
                labSize = 6,
                xlim = c(-5,5),
                ylim = c(-1,5),
                legendPosition = 'right')

3. Plotting top genes on a heatmap:

#Plotting DEGs
ComplexHeatmap::pheatmap(sort.trnsf.df[rownames(sort.trnsf.df) %in% upreg.genes,], scale= "row",
                         cluster_rows = T, cluster_cols = F, annotation_legend = T, show_colnames = F, show_rownames = T,
                         legend = T, treeheight_row = 0, color = heatmap.color.code,
                         annotation_col=annot_df, annotation_colors = annot_colors, fontsize_row = 8)

ComplexHeatmap::pheatmap(sort.trnsf.df[rownames(sort.trnsf.df) %in% downreg.genes,], scale= "row",
                         cluster_rows = T, cluster_cols = F, annotation_legend = T, show_colnames = F, show_rownames = T,
                         legend = T, treeheight_row = 0, color = heatmap.color.code,
                         annotation_col=annot_df, annotation_colors = annot_colors, fontsize_row = 8)

Enrichment testing results

1. List enrichments

#DT::renderDT(data.frame(res.complex), "OH1.5A",scrollY=1000)
paste("Reminder: we have", length(upreg.genes), "upregulated and", length(downreg.genes), "downregulated genes.")

[1] “Reminder: we have 17 upregulated and 8 downregulated genes.”

createDT(DEG.enrich.res, "Enrichment", scrollY=1000)

2. GO enrichments

dotplot(GO.up, split="ONTOLOGY") + facet_grid(ONTOLOGY~., scale="free")
dotplot(GO.down, split="ONTOLOGY") + facet_grid(ONTOLOGY~., scale="free")

Signature co-variance

1. Signature expression (Microglia, NfKB, IFNy)

mic.panel <- mic.panel[,order(int.covs$status, decreasing=F)]
pheatmap::pheatmap(mic.panel, cluster_rows = T, cluster_cols = F, annotation_legend = F, show_colnames = F, show_rownames = T,
                   legend_breaks=c(-4,0,4,4), scale="row",
                   legend = F, treeheight_row = 0,
                   legend_labels=c("-4","0","4","Residual expression \n\n"),
                   annotation_col=annot_df, annotation_colors = annot_colors, fontsize_row = 8)
NFKB.panel <- NFKB.panel[,order(int.covs$status, decreasing=F)]

pheatmap::pheatmap(NFKB.panel, cluster_rows = T, cluster_cols = F, annotation_legend = F, show_colnames = F, show_rownames = F,
                   legend_breaks=c(-4,0,4,4), scale="row",
                   legend = F, treeheight_row = 0,
                   legend_labels=c("-4","0","4","Residual expression \n\n"),
                   annotation_col=annot_df, annotation_colors = annot_colors, fontsize_row = 8)
IFN.panel <- IFN.panel[,order(int.covs$status, decreasing=F)]

pheatmap::pheatmap(IFN.panel, cluster_rows = T, cluster_cols = F, annotation_legend = F, show_colnames = F, show_rownames = F,
                   legend_breaks=c(-4,0,4,4), scale="row",
                   legend = F, treeheight_row = 0,
                   legend_labels=c("-4","0","4","Residual expression \n\n"),
                   annotation_col=annot_df, annotation_colors = annot_colors, fontsize_row = 8)

(1A) Signature expression only for lFC < -1 | lfc > 1 genes (Microglia, NfKB, IFNy)

for (i in c("NFKB", "IFN", "Microglia")){
  panel <- panel.list[[i]]
  DE <- DE.list[[i]]
  df <- panel[DE$log2FoldChange < -0.5 | DE$log2FoldChange > 0.5,]
  df <- df[,order(as.numeric(colnames(df)), method="radix", decreasing=F)]
  annot_df_heat <- data.frame("Status" = int.covs[order(as.numeric(rownames(int.covs), decreasing=F)),]$status, check.names=F)
  rownames(annot_df_heat) <- colnames(df)
  ComplexHeatmap::pheatmap(df, scale= "row",
                         cluster_rows = T, cluster_cols = F, annotation_legend = T, show_colnames = F, show_rownames = T,
                         legend = T, treeheight_row = 0, color = heatmap.color.code,
                         annotation_col=annot_df_heat, fontsize_row = 8)
}

2. Signature co-variance

paste("Microglia")

[1] “Microglia”

PCA_cov_cor_R(int.covs, mic.panel)
PCAplot(mic.pca, "status", PoV.df=mic.PoV, pc.1=2, pc.2=6)

paste("NFKB")

[1] “NFKB”

PCA_cov_cor_R(int.covs, NFKB.panel)
PCAplot(IFN.pca, "status", PoV.df=IFN.PoV, pc.1 = 5, pc.2 = 6)

paste("IFNy")

[1] “IFNy”

PCA_cov_cor_R(int.covs, IFN.panel)
PCAplot(NFKB.pca, "status", PoV.df=NFKB.PoV, pc.1=3, pc.2=6)

(2A) Scatter plot of signature co-variance

scatters.all

(2A) Scatter plot of signature with DEG co-variance (messy, need to clean up)

DEGscatters.all

[[1]]

[[2]]

[[3]]

[[4]]

[[5]]

[[6]]